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ABSTRACT 

The measurement of Doppler velocity shifts in spectra is a ubiquitous theme in astronomy, usuaUy 
handled by computing the cross-correlation of the signals, and finding the location of its maximum. 
This paper addresses the problem of the determination of wavelength or velocity shifts among multiple 
spectra of the same, or very similar, objects. We implement the classical cross-correlation method 
and experiment with several simple models to determine the location of the maximum of the cross- 
correlation function. We propose a new technique, self-improvement, to refine the derived solutions 
by requiring that the relative velocity for any given pair of spectra is consistent with all others. By 
exploiting all available information, spectroscopic surveys involving large numbers of similar objects 
may improve their precision significantly. As an example, we simulate the analysis of a survey of G- 
type stars with the SDSS instrumentation. Applying self-improvement refines relative radial velocities 
by more than 50% at low signal-to-noise ratio. The concept is equally applicable to the problem of 
combining a series of spectroscopic observations of the same object, each with a different Doppler 
velocity or instrument-related offset, into a single spectrum with an enhanced signal-to-noise ratio. 
Subject headings: techniques: radial velocities — methods: statistical — stars: kinematics 



1. INTRODUCTION 

The application of cross-correlation techniques to mea- 
sure velocity shifts has a long history (Simkin 1972, 1974; 
Lacy 1977; Tonry & Davis 1979), and with the advent 
of massive digital spectroscopic surveys of galaxies and 
stars, the subject has renewed interest. The recently 
completed Sloan Digital Sky Survey (SDSS) has col- 
lected spectra for more than 600,000 galaxies and 90,000 
quasars (Adelman-McCarthy et al. 2007, York et al. 
2000). The SDSS has also obtained spectra for about 
200,000 galactic stars, and it is now being extended at 
lower galactic latitudes by SEGUE with at least as many 
spectra (Rockosi 2005, Yanny 2005). Another ongo- 
ing galactic survey, RAVE, is expected to collect high- 
resolution spectra for a million stars by 2011 (Steinmetz 
et al. 2006), and the plans for the GAIA satellite include 
measuring radial velocities for 10^ stars by 2020 (Katz 
et al. 2004). Extracting the maximum possible informa- 
tion from these spectroscopic surveys requires carefully 
designed strategies. 

Cross-correlation has been the target of numerous de- 
velopments in recent years (see, e.g., Mazeh & Zucker 
1994, Statler 1995, Torres, Latham & Stefanik 2007, 
Zucker 2003), but several practical aspects of its imple- 
mentation would benefit from further research. These 
include the selection of templates (e.g., observed vs. syn- 
thetic libraries), how to combine measurements from 
multiple templates, the method to determine the maxi- 
mum of the cross-correlation function, data filtering, and 
error determination. Some of these issues are briefly ad- 
dressed in this paper, but our focus is on how the require- 
ment of coherence among all entries in a radial velocity 
data base can be used to improve the original measure- 
ments. A different but plausible approach has been re- 
cently proposed by Zucker & Mazeh (2006). 

The Doppler shifts of targets in a spectroscopic survey 
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are determined one at a time. Each object's projected ve- 
locity is measured independently, not counting a possible 
common set of cross-correlation templates. For a given 
template, from any pair of (projected) velocity measure- 
ments, we can derive a relative velocity between the two 
objects involved. However, that figure will likely be nu- 
merically different from the value inferred from the direct 
cross-correlation between their spectra, even if the two 
objects are of the same class. In this paper, we argue that 
it is possible to improve the original determinations by 
imposing consistency among all available measurements. 
Our discussion is oriented to the case of a homogeneous 
sample: multiple observations of the same or similar ob- 
jects. 

In the following section I introduce cross-correlation, 
with a brief discussion about error evaluation. Section 
[3] presents the notion of self-improvement and Section 
|4] extends the method to the more realistic scenario in 
which the spectra in a given data set have varying signal- 
to-noise ratios. In f}5] we explore an application of the 
proposed technique involving low-resolution spectra, con- 
cluding the paper with a brief discussion and reffections 
about future work. 

2. CROSS-CORRELATION ANALYSIS 

The most popular procedure for deriving relative ve- 
locities between a stellar spectrum and a template is the 
cross-correlation method (Tonry & Davis 1979). This 
technique makes use of all the available information in 
the two spectra, and has proven to be far superior than 
simply comparing the Doppler shifts between the central 
wavelengths of lines when the signal-to-noise ratio is low. 
The cross-correlation of two arrays (or spectra) T and S 
is defined as a new array C 
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If the spectrum T is identical to S, but shifted by an in- 
teger number of pixels p, the maximum value in the array 
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C will correspond to its element i = p. Cross-correlation 
can be similarly used to measure shifts that correspond 
to non-integer numbers. In this case, finding the location 
of the maximum value of the cross-correlation function 
can be performed with a vast choice of algorithms. 

The most straightforward procedure to estimate real- 
istic uncertainties involves an accurate noise model and 
Monte-Carlo simulations, and that is the method we use 
in Section [5] We employ Gaussians and low-order poly- 
nomials to model the peak of the cross-correlation func- 
tion. For these simple models, implemented in a com- 
panion IDL code, it is possible to derive analytical ap- 
proximations that relate the uncertainty in the location 
of the maximum of the cross-correlation function to the 
covariance matrix [Uy ]. 



will not be reduced that fast. In addition, all Vij are also 
correlated with all Vik, driving the improvement farther 
away from the ideal l/y/n behavior. We can expect that 
after a sufficient number of spectra are included, either 
random errors will shrink below the systematic ones or 
all the available information will already be extracted, 
and no further improvement will be achieved. 

4. GENERAL CASE OF MULTIPLE SPECTRA WITH 
DIFFERENT SIGNAL-TO-NOISE RATIOS 

The case addressed in Section [3] corresponds to a set 
of spectra of the same quality. If the uncertainties in 
the measured relative radial velocities differ significantly 
among pairs of spectra, Eq. [3]can be generalized by using 
a weighted average 



3. VELOCITY SELF-IMPROVEMENT 

Digital cross-correlation, introduced in Section [21 is 
commonly employed to derive Doppler radial velocities 
between two spectra. The discussion in this section is, 
nonetheless, more general, and deals with the statistical 
improvement of a set of relative velocity measurements. 

If three spectra of the same object are available and we 
refer to the relative radial velocity between the first two 
as Vi2, an alternative estimate of V12 can be obtained 
by combining the other relative velocity measurements, 
Vi3 — V23. Assuming uniform uncertainties, the error- 
weighted average of the two values is ¥(2 = {2Vi2 + V13 — 
^23)73. For a set of n spectra, we can obtain an improved 
relative radial velocity determination between the pair ij 
by generalizing this expression 
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It can be seen from Eq. [l]that the correlation of T and 
S is equal to the reverse of the correlation between S and 
T. Thus, when the relative velocities between two spectra 
is derived from cross-correlation and the spectra have a 
common sampling, it will be satisfied that Vij = —Vji, 
but this will not be true in general. For example, if we are 
dealing with grating spectroscopy in air, changes in the 
refraction index with time may alter the wavelength scale 
and the spectral range covered by any particular pixel, 
requiring interpolation. If our choice is to interpolate 
the second spectrum (S) to the scale of the first (T), this 
may introduce a difference between Vij and —Vji due to 
different interpolation errors. We can accommodate the 
general case by writing 
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Note that this definition ensures that V/, 
V- = 0. 



'V'i, and 



If the quality of the spectra is uniform, and all mea- 



sured radial velocities V, 
of the same size a 
an uncertainty a 
different from —Vi 



have independent uncertainties 
aij, the primed values would have 
n. Despite Vij may be numerically 
jj, aij will be highly correlated with 
'ji, and thus the uncertainty in the primed velocities 
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and the uncertainty is 
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In the common case in which Va 



-Vij, the counter- 



part of Eq. [2] for dealing with spectra of varying signal- 
to-noise ratios reduces to 



where 




(7) 



(8) 



In the next section we use simulated spectra for a case 
study: multiple observations of the same object or mas- 
sive surveys involving large numbers of very similar ob- 
jects at intermediate spectral resolution. 

5. EXAMPLE: AN SDSS-LIKE SURVEY OF G-TYPE STARS 

The SDSS spectrographs deliver a resolving power 
of i? = X/FWHM ~ 2000, over the range 381-910 
nm. These two fiber-fed instruments are attached to a 
dedicated 2.5m telescope at Apache Point Observatory 
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(Gunn et al. 2006). Each spectrograph can obtain spec- 
tra for 640 targets simultaneously. As a result of a fixed 
exposure time in SDSS spectroscopic observations, the 
flux in a stellar spectrum at a reference wavelength of 
500 nm, /500, correlates well with the g magnitude of 
the star and with the signal-to-noise ratio at 500 nm 
{S/N^oo)- On average, we find S/N^oo = 25 at 5 ~ 17.85 
mag. To build a realistic noise model, we used the fluxes 
and uncertainties for 10,000 spectra publicly released as 
part of DR2 (Abazajian et al. 2004) to derive, by least- 
squares fitting, a polynomial approximation. When /500 
is expressed in erg cm^ s~^ which are the units used 
in the SDSS data base, we can write 

S/N50Q ~ -167386.14- 41737.014/:-3884.2333/:2 

-160.00277/:^ - 2.4627508/:^ (9) 

where C — logj^Q /500 . This relationship holds in the 
range -16.5 < £ < -14.5. 

The uncertainties in the SDSS fluxes for stars - mostly 
relatively bright calibration sources- are not dominated 
by photon noise, but by a floor noise level of 2-3% asso- 
ciated with a combination of effects, including imperfect 
flat-fielding and scattered light corrections. Errors are 
highly variable with wavelength, but the noise at any 
given wavelength depends linearly on the signal. Based 
on the same set of SDSS spectra used for Eq. [9l we 
determine the coefficients in the relation 

(T\^ax + bxfx, (10) 

which we use here for our numerical experiments. For a 
given choice of S/N500, we interpolate the table of coef- 
ficients ax and bx derived from SDSS data, and by in- 
verting Eq. [9] we derive the fiux at 500 nm. Finally, we 
scale the spectrum fluxes and calculate the expected er- 
rors at all wavelengths using Eq. [101 Gaussian noise is 
introduced for each pixel position, according to the ap- 
propriate error, simulating multiple observations of the 
same star to create an entire library of spectra. 

We employed a spectrum of HD 245, a nearby G2 star^, 
to produce spectra that resemble SDSS observations with 
various signal-to-noise ratios. Radial velocities are also 
artificially introduced. The spectrum of HD 245 used 
here has a resolving power of i? ~ 10, 000^ and is in- 
cluded in the Elodie.3 database (Moultaka et al. 2004, 
Prugniel & Soubiran 2001). As the rest of the library, 
this spectrum was obtained with the 1.9m telescope and 
the Elodie spectrograph at Haute Provence. The origi- 
nal fluxes are resampled to In A, and then smoothed to 
R = 2000 by Gaussian convolution. The output fluxes 
are sampled with 12 pixels per resolution element. The 
Doppler shift due to the actual radial velocity of HD 245 
has already been corrected in the Elodie library. New 
values for the radial velocity in the library of simulated 
SDSS observations are drawn from a Normal distribution 
with a cr = 120 km s~^, as to approximate the typical 
range found in F- and G-type stellar spectra included 
in the SDSS (mostly thick-disk and halo stars). The 

^ Metallicity ([Fe/H]~ -0.7), surface gravity (logg ~ 3.7), and 
kinematics, make this object a prototypical thick-disk turn-ofF star. 

^ The resolving power of this spectrum varies somewhat with 
wavelength from about _R ~ 11000 at 660 nm to 7700 at 480 nm. 
This variation is, however, irrelevant when smoothing the data to 
R = 2000 as we do in these experiments. 



wavelength scale is then 'Doppler' shifted, changed to 
vacuum (Aq), and the spectrum resampled with a step 
of 10"'* in log;^o -^0 (approximately 2.17 pixels per resolu- 
tion element). The Elodie spectra only cover the range 
400 < Aair < 680 nm, and therefore a similar range is 
finally kept for the SDSS-style files, which include 2287 
pixels. 

5.1. Cross- correlation 

We employed a set of 40 test spectra with S/N^qq ~ 50, 
measuring the relative radial velocities for all possible 
pairs. Fig. [1] illustrates two sample spectra and their 
cross-correlation function. To avoid very large or small 
numbers, the input arrays are simply divided by their 
maximum values before cross-correlation. We used sec- 
ond and third order polynomials, as well as a Gaussian 
to model the cross-correlation function and estimate the 
location of its maximum by least-squares fitting. The 
solid line in the lower panels of the figure are the best- 
fitting models. We experimented varying the number of 
pixels involved in the least-squares fittings {Nf). 

With the sampling used, the measured relative shifts 
in pixel space (x) correspond to a velocity V = xc x 
lO^"*/ logj^Q e, where c is the speed of light in vacuum; 
one pixel corresponds to 69 km s~*. We compare the rel- 
ative velocities between all pairs of spectra derived from 
the measurement of the location of the cross-correlation 
peaks with the known, randomly drawn, relative veloc- 
ities. The average difference for the 1600 velocities (40 
spectra) < dV > and the rms scatter (cr) are used to 
quantify systematic and random errors, respectively. 

Our experiments exhibit no systematic errors in the 
derived velocity when the number of points entering the 
fit Nf was an odd number, i.e., when we use the same 
number of data points on each side from the pixel closest 
to the peak of the cross-correlation function. Modest 
offsets (< SV > /cr ~ 0.02), however, are apparent when 
fitting polynomials to an even number of data points, 
despite we enforce the maximum to be bracketed by the 
two central data points. 

Random errors increase sharply with the number of 
data points involved in the fittings for the polynomial 
models, but not for the Gaussian model. The best re- 
sults for the polynomials are obtained when the lowest 
possible orders are used. Using less than 11 points for the 
Gaussian did not produce reliable results, as there was 
not enough information to constrain all the parameters 
of the model, which includes a constant base line. The 
best performance a — 2.8 km s^^ was obtained using a 
second order polynomial and Nf = 3. Using a Gaussian 
model achieved a minimum a = 4.3 km s"*, fairly inde- 
pendent of Nf . The third order polynomial provided the 
poorest performance, a = 5.4 km s~* at best. 

The cross-correlation can be computed in Fourier 
space, taking advantage of the correlation theorem 
(Brigham 1974). This fact is usually exploited to speed 
up the calculation dramatically, as fast Fourier trans- 
forms can be calculated with a number of operations pro- 
portional to A^log2 N, compared to iV^ required by Eq. 
[TJ Note, however, that for medium-resolution surveys of 
galactic stars, the velocity offsets, limited by the galactic 
escape velocity, usually correspond to a limited number 
of pixels. Therefore, it is only necessary to compute the 
values of C in the vicinity of the center of the array, ren- 
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Fig. 1. — Two sample spectra used in the tests described in i|5.1l They correspond to the same original spectrum, shifted by an arbitrary 
radial velocity, and with added Gaussian noise so that S/N^qq = 50. (The fluxes have been also arbitrarily scaled for display). The inset 
in the upper panel is a close-up of the blue part of the spectrum. The lower panels show the cross-correlation vector and illustrate the 
determination of the maximum by using Gaussian, quadratic, and cubic models with Nf = 20, 7, and 8 data points, respectively. The 
dashed vertical lines show the location of the maxima, and the horizontal dashed lines indicate the estimated uncertainties. 

dering the timing for a direct calculation similar to one 
performed in transformed space'^. 

5.2. Self-improvement 

To test the potential of the proposed self-improvement 
technique we repeat the same exercise described in 
but using increasingly larger datasets including up to 320 
spectra, and adopting three different values for the S/N 
per pixel at 500 nm: 50, 25, and 12.5. We calculated 
the cross-correlation between all pairs of spectra (matrix 
V), and performed quadratic fittings to the 3 central data 
points, cubic polynomial fittings to the central 4 points, 
and Gaussian fittings involving the 11 central points. We 
estimated the uncertainties in our measurements by cal- 
culating the rms scatter between the derived and the 
known relative velocities for all pairs. Then we applied 
Eq. [2] to produce a second set of self-improved veloci- 

In our tests with a Sun Fire V210 (1 CPU) and IDL 6.1, the 
timings were similar in pixel or Fourier space: a a single cross- 
correlation (N = 2287) took 2.5 — 2.6 X 10"'^ seconds in Fourier 
space (arrays padded to 5^ or 2^^), while in pixel space, with a lag 
range restricted to ±15 pixels (±1035 km s~^), it took 2.2 X 10~^ 
seconds. 




Fig. 2. — Distribution of errors in the relative radial velocities 
derived from cross-correlation (quadratic model) for the case of 
n = 120 spectra and S/N^oo = 25. The solid line represents the 
original distribution, and the dashed line the result after applying 
self-improvement. The error distributions are symmetric because 
the array V is antisymmetric. 
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Fig. 3. — Estimated Icr uncertainties in the cross-correlation between all pairs in a set of noise-injected spectra of HD 245 as a function of 
the number of spectra in) for three values of the signal-to-noise ratio at 500 nm. Three models are considered to determine the maximum 
of the cross-correlation function. The black lines correspond to the original measurements, and the red lines to the final uncertainties after 
self-improvement. 

second random number generator'^ and the results were 
consistent. 

As described in ijS.ll the quadratic model per- 
forms better on the original velocity measurements for 
5'/^50o = 50 and 25. At the lowest considered S/N value 
of 12.5, however, the Gaussian model delivers more accu- 
rate measurements. Self-improvement reduces the errors 
in all cases. Although a second order polynomial fitting 
works better than third order for the original measure- 
ments, the two models deliver a similar performance af- 
ter self-improvement. Interestingly, the impact of self im- 
provement is smaller on the results from Gaussian fittings 
than on those from polynomial fittings. As expected, the 
errors in the original measurements are nearly indepen- 
dent of the number of spectra in the test, but there is 
indication that at low signal-to-noisc the errors after self- 
improvement for the polynomial models decrease as the 
sample increases in size, until they plateau for n > 100. 

From these experiments, we estimate that the best ac- 
curacy attainable with the original cross-correlation mea- 
surements are about 3, 6, and 15 km s~-^ at S/N^oo — 50, 
25, and 12.5, respectively. Our results also indicate that 
by applying self-improvement to samples of a few hun- 
dred spectra, these figures could improve to roughly 2.5, 



common to all spectra, the matrix V is antisymmetric 
and we can use Eq. [2] instead of Eq. [3l) A first effect of 
the transformation from V to V, is that the systematic 
offsets described in !j5.1l when using polynomial fittings 
with even values of TV/ disappear (the same systematic 
error takes place for measuring Vij and Vji, canceling out 
when computing V^'j). 

More interesting are the effects on the width of the 
error distributions. Fig. [2] illustrates the case when a 
quadratic model is used for n ~ 120 and S/N500 — 25. 
The solid line represents the original error distribution 
and the dashed line the resulting distribution after self- 
improvement. Fig. [3] shows the rms scatter as a func- 
tion of the number of spectra for our three values of 
the S/N ratio at 500 nm. The black lines lines show 
the original results, and the red lines those obtained af- 
ter self-improvement. Each panel shows three sets of 
lines: solid for the quadratic model, dotted for the cu- 
bic, and dashed for the Gaussian. Extreme outliers at 
\V\ > 2000(5/^^500)"^ km s"^, if any, were removed be- 
fore computing the width of the error distribution (cr). 
Note the change in the vertical scale for the case with 
S/N5Q0 — 12.5. For the experiments with n < 80, several 
runs were performed in order to improve the statistics, 
and the uncertainty (standard error of the mean) is indi- 
cated by the error bars. These results are based on the 
Gaussian random-number generation routine included in 
IDL 6.1, but all the experiments were repeated with a 



■* The F90 module random. f 90 by Alan J. Miller, and in par- 
ticular the normal random number algorithm presented by Leva 
(1992), coupled to the intrinsic raiidom_number function included 
with the Sun Workshop 2.0 F90 compiler. 
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4, and 9 km s^^ at S/N500 — 50, 25, and 12.5, respec- 
tively. 

We obtained an independent estimate of the precision 
achievable by simply measuring for 320 spectra the wave- 
length shift of the core of several strong lines (Ha, H/3, 
H7, and US) relative to those measured in the solar spec- 
trum (see Allende Prieto et al. 2006), concluding that 
radial velocities can be determined from line wavelength 
shifts with a la uncertainty of 3.8 km s~^ at S/N = 50, 
7.2 km s"^ at S/N = 25, and 15.9 km s'^ at S/N = 12.5 
- only 10-20% worse than straight cross-correlation- but 
these absolute measurements cannot take advantage of 
the self-improvement technique. 

Allende Prieto et al. (2006) compared radial veloci- 
ties determined from the wavelength shifts of strong lines 
for SDSS DR3 spectra of G and F-type dwarfs with the 
SDSS pipeline measurements based on cross-correlation. 
The derived la scatter between the two methods was 
12 km or, assuming similar performances, a preci- 
sion of 8.5 km s~^ for a given method^. The spectra 
employed in their analysis have a S/N^qq distribution 
approximately linear between S/N^qq = 10 and 65, with 
dN / d{S / N) ^ 7.7 and with mean and median values of 
22 and 18, respectively. Their result is in line with the 
expectations based on our numerical tests that indicate a 
potential precision of 6-7 km s^^ at S/N^qq = 25. Inde- 
pendent estimates by the SDSS team are also consistent 
with these values (Rockosi 2006; see also www . sdss . org). 

After correcting for effects such as telescope flexure, 
the wavelength scale for stellar spectra in DR5 is accu- 
rate to better than 5 km s~^ (Adelman-McCarthy et al. 
2007). This value, derived from the analysis of repeated 
observations for a set of standards and from bright stars 
in the old open cluster M67, sets an upper limit to the 
accuracy of the radial velocities from SDSS spectra, but 
random errors prevail for S/N < 25. Provided no other 
source of systematic errors is present, our tests indicate 
that self-improvement could reduce substantially the typ- 
ical error bars of radial velocities from low signal-to- noise 
SDSS observations. 

6. DISCUSSION AND CONCLUSIONS 

This paper deals with the measurement of relative 
Doppler shifts among a set of spectra of the same or 
similar objects. If random errors limit the accuracy of 
the measured relative velocity between any two spec- 
tra, there is potential for improvement by enforcing self- 
consistency among all possible pairs. 

This situation arises, for example, when a set of spec- 
troscopic observations of the same object are available 
and we wish to co-add them to increase the signal-to- 
noise ratio. The spectra may be offset due to Doppler 
velocity offsets or instrumental effects, the only differ- 
ence being that in the former case the spectra should be 
sampled uniformly in velocity (or log A) space for cross- 
correlation, while in the latter a different axis may be 
more appropriate. 

Another application emerges in the context of surveys 
that involve significant numbers of spectra of similar ob- 

^ Note that this figure was determined by fitting a Gaussian 
profile to the distribution of difli'erences, and it will increase some- 
what if a straight rms difference is considered instead. On the other 
hand, the spectral types of the stars considered spanned a broader 
range (F and G) than in our simulations. 



jects. Radial velocities for individual objects can be 
derived using a small set of templates and later self- 
improved by determining the relative velocities among all 
the survey targets and requiring consistency among all 
measurements. The potential of this technique is illus- 
trated by simulating spectra for a fictitious survey of G- 
type turn-off stars with the SDSS instrumentation. Our 
simulations show that applying self-improvement has a 
significant impact on the potential accuracy of the deter- 
mined radial velocities. The tests performed dealt with 
relative velocities, but once the measurements are linked 
to an absolute scale by introducing a set of well-known 
radial velocity standards in the sample, the relative val- 
ues directly translate into absolute measurements. The 
ongoing SEGUE survey includes, in fact, large numbers 
of G-type stars, and therefore our results have practical 
implications for this project. 

The proposed scheme handles naturally the case when 
multiple templates are available. Templates and targets 
are not treated differently. Relative velocities are mea- 
sured for each possible pair to build Vij , and consistency 
is imposed to derive V/j by using Eqs. [3] or ID If, for ex- 
ample, the templates have been corrected for their own 
velocities and are the first 10 spectra in the sample, the 
velocity for the j— th star {j > 10) can be readily ob- 
tained as the weighted average of the V/j elements, where 
i runs from 1 to 10. The final velocities would take advan- 
tage of all the available spectra, not just the radial veloc- 
ity templates, with differences in signal-to-noise among 
spectra already accounted for automatically. 

Very recently, Zucker & Mazeh (2006) have proposed 
another approach with the same goals as the method 
discussed here. Their procedure determines the rela- 
tive velocities of a set of n spectra by searching for the 
Doppler shifts that maximize the value of the parame- 
ter p = (Am — 1)/('T- — l)j where Xm is the maximum 
eigenvalue of the correlation matrix - a two-dimensional 
array whose ij element is is the cross-correlation function 
between spectra i and j. Zucker & Mazeh's algorithm is 
quite different from the self-improvement method pre- 
sented here. It involves finding the set of velocities 
that optimally aligns the sample spectra, whereas self- 
improvement consists on performing very simple alge- 
braic operations on a set of radial velocities that have 
already been measured. Self-improvement is obviously 
more simple to implement, but a detailed comparison 
between the performance of the two algorithms in prac- 
tical situations would be very interesting. 

This paper also touches on the issue of error deter- 
mination for relative radial velocities derived from cross- 
correlation, and convenient analytical expressions are im- 
plemented in an IDL code available online^. We have 
not addressed many other elements that can potentially 
impact the accuracy of Doppler velocities from cross- 
correlation, such as systematic errors, filtering, sam- 
pling, or template selection. The vast number of spectra 
collected by current and planned spectroscopic surveys 
should stimulate further thought on these and other is- 
sues with the goal of improving radial velocity determi- 
nations. There is certainly an abundance of choices that 
need to be made wisely. 

^ To gether with the electronic version of this article or from 
|http://hebe.as. utexas.edu/stoois7l 
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